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ABSTRACT 

In dilute astrophysical plasmas, thermal conduction is primarily along magnetic field lines, and 
therefore highly anisotropic. As a result, the usual convective stability criterion is modified from a 
condition on entropy to a condition on temperature. For small magnetic fields or small wavenumbers, 
instability occurs in any atmosphere where the temperature and pressure gradients point in the same 
direction. We refer to the resulting convective instability as the magnetothermal instability (MTI). 
We present fully three-dimensional simulations of the MTI and show that saturation results in an 
atmosphere with different vertical structure, dependent upon the boundary conditions. When the 
temperature at the boundary of the unstable layer is allowed to vary, the temperature gradient relaxes 
until the unstable region is almost isothermal. When the temperature at the boundary of the unstable 
region is fixed, the magnetic field is reoriented to an almost vertical geometry as a result of buoyant 
motions. This case exhibits more vigorous turbulence. In both cases the resulting saturated heat flux 
is almost one-half of the value expected if the conduction were purely isotropic, Q ~ x^T/L, where 
X is the thermal conductivity, AT is the fixed temperature drop across the simulation domain, and L 
is the temperature gradient scale length. The action of the MTI results in dynamical processes that 
lead to significant transport perpendicular to the initial direction of the magnetic field. The resulting 
magnetoconvection in both cases amplifies the magnetic field until it is almost in equipartition with 
sustained subsonic turbulence. These results are relevant to understanding measurements of the 
temperature profiles of the intracluster medium of clusters of galaxies as well as the structure of 
radiativcly inefficient accretion flows. 

Subject headings: accretion, accretion disks — convection — instabilities — MHD — galaxies: clusters 
— stars: neutron 
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1. INTRODUCTION 

In many astrophysical plasmas the collision frequency 
is much smaller than the Larmor frequency, so that par- 
ticles spiral around the magnetic field lines for very long 
distances between collisions. As a result, transport of 
heat and momentum by thermal conduction and viscos- 
ity is highly anisotropic with respect to the magnetic 
field orientation. In a dilute magnetized plasma, elec- 
tron thermal conductivity parallel to the magnetic field 
is many orders of magnitude larger than both the per- 
pendicular electron conductivity and either component 
of the ion conductivity provided u> < v ee <C £^ e , where 
v ee is the electron collision frequency, Q e is the electron 
gyrofrequency, and u> is the frequency of a physical pro- 
cess of interest. In this regime the plasma is decribed 
by the equations of magnetohydrodynamics (MHD) with 
the addition of B raginskii anisotropic transport terms 
(|Braginskii Ill965t ). A Braginskii description is valid for 
describing processes where the mean free path is a sub- 
stantial fraction of the size of the system, but does not 
permit an analysis of processes that occur on lengthscales 
shorter than the mean free path, e.g. Landau damping. 

These properties can have profound effects on the 
resulting plasma dynamics. A striki ng example is 
the magnetothermal i nstability (MTI) (jBalbus I l2000t 
IParrish fe Stonel l2005h . The magnetothermal instabil- 
ity occurs as a result of heat flow parallel to the mag- 
netic field in a statified atmosphere. The criterion 

1 Program in Applied and Computational Mathematics, Prince- 
ton University, Princeton, NJ, 08544 



for convective stability is modified from the well-known 
Schwarzschild criterion (dS/dz > 0), where the entropy, 
S = VP 1 > to the Balbus criterion (dhiT/dz > 0), where 
p and T are the pressure and temperature respectively. 
Atmospheres that have temperature (as opposed to en- 
tropy) profiles decreasing upwards are unstable to the 
MTI. In this situation the MTI is able to use the tem- 
perature gradient as a source of free energy for the insta- 
bility. The growth rate and saturation of this instability 
in two dimensions has been explored in IParrish fe Stonel 
(|2005h . hereafter PS. In PS we demonstrated that com- 
putationally measured growth rates match theoretical 
estimates from a WKB theory for a variety of bound- 
ary and initial conditions. In addition, we showed that 
MTI-unstable plasmas can produce vigorous convective 
motions and a significant advective heat flux. The satu- 
rated state for an adiabatic atmosphere in 2D was shown 
to result in an isothermal atmosphere. 

We are motivated to study the nonlinear regime of the 
magnetothermal instability in three dimensions for sev- 
eral reasons. First, the saturation mechanism in three 
dimensions may be different from two dimensions. Sec- 
ond, there is a potential for a magnetic dynamo in three 
dimensions. Finally, the nature of convection is known 
to be intrinsically different in three dimensions. 

In this paper we present the saturation and heat trans- 
port properties of the MTI in three dimensions. We 
find the instability is capable of generating a magnetic 
dynamo that greatly amplifies the initial field until it 
is roughly in equipartition with the fluid turbulent ki- 
netic energy. To estimate the efficiency of heat trans- 
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port, we can compare to a fiducial heat flux that would 
be expected if the conduction were purely isotropic, 
Q ~ x^T/L, where x is the thermal conductivity, AT is 
the fixed temperature drop across the simulation domain, 
and L is the temperature gradient scale length. We find 
that the saturated state of the atmosphere varies depend- 
ing on the boundary condition, but is consistently able 
to transport a large fraction (30-50%) of the fiducial heat 
flux. The understanding of the three dimensional behav- 
ior of the MTI is likely to be very applicable to several 
astrophysical systems, such as the intracluster medium 
of clusters of galaxies, the atmospheres of neutron stars, 
and radiatively inefficient accretion flows. 

This paper is organized as follows. In §S]we describe 
our numerical method and the physical mechanism re- 
sponsible for the MTI. In §S] we briefly comment on 
the linear properties of the MTI before focusing on the 
three-dimensional saturation mechanism in three differ- 
ent principal cases. In §3] we discuss several potential 
applications of the MTI. Finally, in Sj5]we summarize our 
results and discuss future work. 

2. METHOD AND PHYSICS OF THE MTI 

2.1. Equations of MHD with Anisotropic Heat 
Conduction 

We solve the usual equations of ideal MHD with the 
addition of a vector heat flux, Q, and a vertical gravita- 
tional acceleration g — ~go&, 
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where the symbols have their usual meaning. The total 
energy E is given as 



E = e- 



B B 



(5) 



and the internal energy, e = p/(j — 1)- For this paper, 
we assume 7 = 5/3 throughout. Note that equations 
([I])-© are in conservative form with the exception of 
the gravitational source terms. The induction equation 
(Eqn. [4]) requires special treatment in MHD to preserve 
the V • B = constraint exactly. 

The heat flux, Q = Qc + Qr, where Qc, the 
anisotropic term is due to Coulombic conduction par- 
allel to magnetic field lines. Cross-field collisions could 
contribute to the isotropic term, but will be neglected 
in this paper due to their very small magnitude in the 
astrophysical situations of interest-the result of the gy- 
roradius being much smaller than the mean free path. 
The isotropic heat flux, Qr is due to radiative trans- 
port. Then 

Q c = -xebb ■ VT, (6) 



Qi 



(7) 



where xc is the Spitzer Coulombic conductivity (|Spitzer I 
1962), b is a unit vector in the direction of the magnetic 
field, and XR is the coefficient of isotropic conductivity. 
For the purposes of the general study of the MTI, we 
consider both coefficients as free parameters. 

2.2. Physics of the MTI 

The physics of the MTI can be understood simply by 
considering a small fluid parcel in a dilute, equilibrium, 
convectively stable at mosphere where the tem perature 
decreases with height (jParrish fe Stone I l2006h . Imag- 
ine a parcel of fluid at height z\ in an equilibrium state 
with the initial magnetic field perpendicular to the tem- 
perature gradient, say in the x direction. In this initial 
configuration b ■ VT = 0, and no heat flows. 

Now, perturb the parcel of fluid to a height z 2 — zi+Sz. 
The displaced parcel adiabatically cools, but the fluid 
element is still connected by the magnetic field to its 
hotter, original depth. The resultant upward heat trans- 
port along the field causes the displaced parcel to become 
buoyant. To understand the instability criterion, we note 
that this discussion requires mechanical equilibrium (me- 
diated by s ound waves) to operate fast er than thermal 
conduction (IChandran fc Dennis II2006T ). The tempera- 
ture of the displaced fluid parcel is that of the lower fluid 
el ement, but its pressur e P* = Pi- In a similar manner 
to IChandran fc Dennis I (|2006h . if we equate these pres- 
sures, assuming P = pT, and utilize a Taylor series for 
T2 about Zi, wc find 



(p* ~ P2)T\ = Szp 2 



dT 
dP 



(8) 



where p* is the density of the displaced element. If the 
density contrast, (p* — p 2 ) < 0, then the fluid parcel 
is buoyant. This realization implies that the instability 
criterion in its simplest form is simply 



dT 

-j— < (instability). 



(9) 



We can be more precise about the instability crite- 
rion and growth rates by following a WKB analysis per- 
formed in detail in § 4. of iBalbus I (|2000l ) and §§3-4 of 
IChandran fc Dennis I (|2006T K We introduce the Brunt- 
Vaisala frequency, N, 

do) 

7/9 oz oz 

the frequency of vertical oscillations in a stably stratified 
atmosphere. Preserving Balbus' notation, we also define 
two useful quantities, 



x'c = —j5-Xc and x = — ^~ {Xc + Xr) ■ (H) 

Using the Fourier convention for perturbations, 
exp {at + ikx), and in the limit of a weak magnetic 
field, the dispersion relation simplifies to the non- 
dimensionalized form 

(£.\ 3 +- (-Y ( x ' Tk2 \ ( a \ dlnr ( x - Tk2 

\NJ 7 VA/7 V N J \nJ dhxS\ N 

(12) 

By analysis of the Routh-Hurwitz criterion, the general 
instability criterion is found to be 
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where the second term is the isothermal limit (7 — ► 1) of 
the Brunt- Vaisala frequency. Examination reveals that 
long wavelengths are the most easily destabilized and 
that strong magnetic fields can suppress the instability. 
A large isotropic conductivity also can short-circuit the 
parallel heat conduction and suppress the instability. In 
the limit of infinite wavelength or very weak magnetic 
field, the instability criterion reduces to the heuristically 
derived result in equation In this limit, any atmo- 
sphere where temperature decreases with height is con- 
vectively unstable to the MTI. 

2.3. Initial Conditions 

For testing linear growth rates, we use the same initial 
vertical equilibrium state outlined in PS §§2.2-2.4. For 
studying the nonlinear saturation properties of the insta- 
bility, we construct an atmosphere similar to §5.2 of PS 
where an MTI-unstable layer is surrounded on top and 
bottom by convectively stable layers. This construction 
prevents the formation of unresolved thermal boundary 
layers at the edge of the computational domain. In par- 
ticular, we begin with the Ansatz that the atmosphere 
exists in a constant gravitational field, g(z) = —goZ, in 
a plane parallel atmosphere. We divide the vertical size, 
L z , into three equal regions. The bottom and top layers, 
region 1 (0 < z < L z /3) and region 3 (2L z /3 < z < L z ) 
respectively, are isothermal atmospheres. In such a re- 
gion, the pressure and density vary exponentially with 
the height. For a region with temperature To, the atmo- 
sphere has solutions of 

p(z) = Po e-^' T \ 



P{z) = P Q c 



-g z/T 



(14) 
(15) 

where Pq = TqPq. Simultaneously, we set the conductiv- 
ity in this region to be purely isotropic with xr = Xiso,o- 
As a consequence of this choice, the second term of equa- 
tion (fT3|) is zero and the MTI is stabilized. For simplicity, 
we choose the parameters at the bottom of the isother- 
mal region to have Tb = Pb = Pb = 1- We designate the 
parameters at top of region 1 with the subscript 1, e.g. 
Tt. 

In the middle region, which we refer to as region 2, 
we construct an atmosphere where the temperature de- 
creases linearly with height as in the previous 2D study. 
Thus our initial Ansatz is a power law atmosphere in the 
region L z /3 < z < 2L z /3: 



T{z) = Ti 



p(z) = Pl 



p(z) =pi 



1 



1 - 



1 - 



(z - L z /3) 



(z - L z /3) 



(z - L z /3) 



(16) 



(17) 



(18) 



where zo is a constant. For stability to adiabatic con- 
vection, we require that < Zq 1 < 2/5go- Our choice 
of zo = 3 satisfies this constraint. Region 2 has conduc- 
tivity that is exclusively anisotropic with xc = Xaniso.o- 
Thus region 2 is unstable to the MTI. We designate the 
hydrodynamic parameters at the top of this region with 
the subscript 2, e.g. T 2 . 



The top layer of the atmosphere, region 3, is con- 
structed in an exactly analogous manner to Region 1. 
It is an isothermal atmosphere with constant tempera- 
ture T t = T3. The thermal conductivity is again set to 
be purely isotropic. 

The magnetic field is initialized only in Region 2. By 
initializing the magnetic field only in the MTI-unstable 
region, its propagation can be used as a tracer of phe- 
nomena like penetrative convection and magnetic pump- 
ing into the stable regions. To initialize the field with a 
fully three dimensional structure, we choose in the region 
L z /3 < z < 2L z /3, 

.(z-LJ3)' 



B(Z) = 




(19) 



y- 



where v = 3 is a half rotation of the helix and represents a 
net flux condition and v = 6 is a full rotation of the helix 
and represents a no net flux condition. The coefficient 
B = 10~ 5 is chosen so that magnetic tension effects are 
unimportant in the linear regime. 

2.4. Dimensionless Analysis 

There are four characteristic time-scales of the prob- 
lem. The first is the sound crossing time of the domain 

t s = L/c B , cl^ 1 —. (20) 
P 

For our standard parameters, the adiabatic soundspeed, 
c s w 1.28, and the domain size is, L = 0.1. Thus, the 
sound crossing time is a short t s m 7.8 x 10 -2 . Likewise, 
we can define an Alfven crossing time as 



T A = L/v A , 



Bi 
4irp 



(21) 



For a typical initial magnetic field of Bo = 10~ 5 , ta ~ 
10 4 , many orders of magnitude longer than the sound- 
crossing time. Thus, we see that Alfven waves are or- 
dered out of the problem, as one would expect for a 
high-/3 plasma. The third timescale is that of the at- 
mospheric oscillation time. This time-scale is simply the 
inverse of the Brunt- Vaisala frequency, 

tn = N~\ (22) 
where N is defined in equation (JTUJ) . For the gradients 
in region 2, N 0.26, which gives tn ~ 3.9. Finally, 
the last timescale of interest is the heat diffusion time 
across the simulation domain. It is obtained by non- 
dimcnsionalizing the thermal conduction equation, re- 
sulting in a time-scale 



L 2 AT 

XcT ' 



(23) 



where AT is the temperature difference across the un- 
stable region, and T is the average temperature. For the 
anisotropic conductivity of xc — 10~ 4 and the standard 
atmosphere, the heat conduction timescale, t x w 1.0. 
The asymptotic ordering of these terms is instructive: 

Ts < t x < t n < ta- (24) 
First, the sound crossing time is rapid compared to the 
buoyancy time. Second, we find that thermal conduc- 
tion is able to smooth the temperature differences on a 
timescale comparable to the buoyant time. Finally, as 
previously noted, the magnetic field is not dynamically 
important. 
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2.5. Numerical Methods 

We use the 3D version of the Athena MHD code 
(jGardiner fc Stone Il2005l |2006[ ) to simulate the MTI far 
into the nonlinear regime. Athena is a higher order, 
fully unsplit Godunov method that combines corner tran- 
port upwind (CTU) methods with constrained transport 
(CT) to evolve the equations of MHD while manifestly 
maintaining the divergence-free constraint of the mag- 
netic field. We utilize the 6 step CTU+CT algorithm. 
Isotropic and anisotropic heat fluxes are calculated as 
face-centered, area-averaged quantities as described in 
detail in the appendix of PS. By finite-differencing the 
anisotropic diffusion terms we are able to keep the ratio of 
parallel conductivity to numerical perpendicular conduc- 
tivity to better than 10 3 for an arbitrary field orientation. 
This level of truncation error is more than sufficient to 
simulate the MTI. Sharma and Hammett (2006, personal 
communication) have found that for certain discontinous 
test problems this differencing can lead to negative tem- 
peratures; however, for our applications the temperature 
gradients are continuous and we observe post hoc that 
no such problems arise. The heat transport algorithm 
is operator split from the MHD and sub-cycled when 
the heat conduction timestep is shorter than the MHD 
Courant condition. Both the MHD and heat conduction 
algorithms are parallelized using MPI. 

All of the non-linear layered runs in this paper are 
performed on a three-dimensional grid initialized as de- 
scribed in 32J3] with gravity in the z direction. The size 



of the domain is L x = L y = 0.1 and L z = 0.2. Our 
standard numerical resolution is 64 2 x 128, such that 
our zones are cubic. We have performed some simula- 
tions at twice this resolution in each direction. Our box 
size is small compared to the atmospheric scale height, 
H ~ c 2 s /g « 1.64, where c s « 1.28 is the adiabatic sound 
speed. Since the MTI is a local instability, the evolution 
of modes with wavelengths much smaller than the scale 
height are independent of the exact form of g(z). 

We use periodic boundary conditions for the MHD 
variables and the temperature in the x and y directions. 
In the vertical direction, we use a modified reflecting 
boundary condition in which all MHD variables are re- 
flected except the internal energy. Due to the non-zero 
gravitational source term at the boundary, we extrapo- 
late the pressure across the boundary to ensure the pres- 
sure gradient balances gravity, and the boundary cell face 
does not undergo spurious acceleration. For the temper- 
ature in the z-direction wc implement a Dirichlct bound- 
ary condition; that is, the temperature is held fixed at the 
boundary and heat is permitted to flow as necessary. We 
refer to this boundary condition as "conducting bound- 
aries." 

To seed multiple modes of the MTI, we add small 
Gaussian white noise perturbations in all three veloc- 
ity components in the unstable region. The standard 
perturbation amplitude is vq = 10~ 4 , which is highly 
subsonic. The standard initial magnetic field strength is 
-B /(47r) 1 / 2 = 10~ D . The Alfven velocity is then roughly 
va ~ 10 -5 (varying slightly with density), which is sim- 
ilar in magnitude to the perturbation velocity, although 
this fact does not affect the instability, as was shown in 



3. SATURATION OF THE MTI IN 3D 

We have tested our 3D numerical algorithms by con- 
firming the linear growth rates of the MTI. As expected, 
we find that the measured growth rates are very close to 
the ones predicted from the WKB theory, and no devi- 
ation is found in the saturation of single-mode pertur- 
bations from the two-dimensional results. We therefore 
quickly turn our attention to the non-linear regime in 
3D. 

Table Q] lists the runs that will be discussed in this 
paper. All of them have the computational set-up de- 
scribed in ij2.3l with constant conductivities. Runs Fl - 
F3 will be our fiducial cases of interest and will be dis- 
cussed in the following sections. These runs have a fixed 
anisotropic thermal conductivity, Xaniso,o = 1 x 10~ 4 , in 
the unstable layer with varying isotropic conductivities 
in the stable layers. Runs Nl - N2 examine the effect 
of a net magnetic flux on the fiducial runs, and run N3 
is a high resolution version of run F3. Runs N4 - N6 
have a higher magnitude of conductivity and are inter- 
esting for comparison to the canonical cases. It should 
be noted that not all of these latter runs have reached 
full saturation, therefore we report the parameters of the 
fiducial runs at an early time for comparison purposes. 
We realized that the saturation mechanism had two very 
different forms with a continuum of intermediate cases. 
The variation in saturation mechanism is directly related 
to the ratio of the conductivity in the MTI-stable outer 
layers (xi so ,o) to the conductivity in the central MTI- 
unstable layer (x an iso,o)- This condition is essentially a 
ratio of the applied heat flux across the unstable layer to 
the parallel conductivity. 



3.1. Case A: Temperature Gradient Relaxation 

We begin our discussion of saturation by examing run 
Fl, in which the conductivity in the isotropic, stable lay- 
ers is an order of magnitude smaller than the anisotropic 
conductivity. Figure [T] shows the evolution of the volume- 
averaged components of the magnetic energy density and 
kinetic energy density over the entire domain, including 
the stable regions. There is a very rapid spike in all 
the component of both magnetic and kinetic energy as 
the linear phase of the instability progresses. The ki- 
netic energy then decays to a relatively constant value 
as the magnetic field is amplified further. Figure [2] 
shows the increase of the volume-averaged Alfven ve- 
locity, (Va) = {{B 2 1 p) 1 ^ 2 ^ , as the magnetic field in- 
creases to (Va) = 4.9 x 10~ 5 . Meanwhile the RMS Mach 
number, (v 2 /c 2 ) , has saturated at a relatively con- 
stant 5.7 x 10~ 5 which corresponds to a vertical velocity, 
(v z ) = (c s )(Ma) w 7.3 x 10~ 5 . This Mach number rep- 
resents highly subsonic turbulence. Thus, the magnetic 
energy density has been amplified such that it is almost 
in equipartition with the kinetic energy density. We can 
quantitatively estimate the increase in magnetic energy 
density as listed in Table [21 with two similar but slightly 
different measures of the magnetic energy based on the 
square of the mean field or the mean of the squared field: 
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TABLE 1 

Properties of Nonlinear Runs 



Run Xiso,0 Xaniso,0 

B /(4tt) 1/2 Flux a Resolution 



Fl 1 x 1CT 5 1 x KT 4 1 x 1CT 5 64 x 64 x 128 

F2 1 x 10" 4 1 x 10" 4 1 x 1CT 5 64 x 64 x 128 

F3 1 x 10~ 3 1 x 10" 4 1 x 1CT 5 64 x 64 x 128 

Nl 1 x 10~ 5 1 x 10" 4 1 x 1CT 5 + 64 x 64 x 128 

N2 1 x 10~ 3 1 x 10~ 4 1 x 1CT 5 + 64 x 64 x 128 

N3 i x io~ 3 i x icr 4 i x icr 5 o 128 x 128 x 256 

N4 5 x 10~ 5 5 x Kr 4 1 x 10" 5 64 x 64 x 128 

N5 5 x 10- 6 5 x 10" 4 1 x 10~ s + 64 x 64 x 128 

N6 5 x 1CT 4 5 x 10" 5 1 x 1CT 5 + 64 x 64 x 128 



a '+' indicates net magnetic flux threading the domain; '0' indicates zero net magnetic flux. 



TABLE 2 

Saturation Properties in the Nonlinear Regime 



Run 


Case 




8(Bf 


S(B 2 ) 


RMS Mach 


(Va) 


(Ob) 




(Q)/Q 


Fl 


A 


2,050 


20.6 


33.4 


8.9 x 10 -5 


4.6 x 10~ 5 


17.5° 


15.8% 


27% 






4, 430* 


22.7 


36.2 


5.7 x 10" 5 


4.9 x 10~ 5 


18.9° 


10.8% 


49% 


F2 


hit 


2,050 


26.7 


45.1 


6.6 x 10~ 4 


5.3 x 10~ 5 


19.8° 


46.9% 


27% 






10, 700* 


186 


298 


6.1 x 10" 4 


1.4 x 10~ 4 


36.6° 


46.0% 


27% 


F3 


B 


2,050 


71.6 


121 


1.8 x 10~ 3 


8.6 x 10~ 5 


25.3° 


90.7% 


45% 






8, 470* 


5290 


7760 


1.5 x 10 -3 


7.5 x 10" 4 


54.9° 


89.1% 


50% 


Nl 


A 


2,050 


74.2 


116 


2.2 x 10~ 4 


8.8 x 10~ 5 


27.1° 


14% 


30% 


N2 


B 


2,050 


379 


615 


1.9 x 10" 3 


2.0 x 10~ 4 


38.7° 


90.7% 


45% 


N3 


B 


2,050 


1680 


3480 


6.1 x 10" 4 


4.1 x 10~ 4 


35.7° 


90.9% 


42% 


N4 


A 


2,050 


22.1 


36.8 


3.5 x 10" 4 


4.8 x 10~ 5 


17.9° 


10.7% 


44% 


N5 


A 


2,050 


97.8 


163 


3.5 x 10~ 4 


1.01 x 10~ 4 


29.5° 


10.6% 


45% 


N6 


B 


2,050 


144 


213 


1.2 x 10~ 3 


1.2 x 10~ 4 


32.2° 


89.5% 


41% 



( a ) The symbol * indicates that the run is considered to be approximately saturated. ( b ) The temperature gradient is given as a percentage 
of the initial temperature gradient in the central region. 



At saturation, run Fl achieves a modest increase in these 
two measures of a factor of 20.6 and 33.4, respectively. 

The evolution to saturation can be understood by ex- 
amining the evolution of the horizontally-averaged tem- 
perature profile shown in Figure [3] The indicated times 
are shown non-dimensionalized to the sound crossing 
time. The temperature profile evolves very rapidly to- 
wards an essentially isothermal atmosphere in the MTF 
unstable, central region. The atmosphere in the stable 
layer slowly relaxes to an essentially linear temperature 
profile. In this case, the temperature profile has relaxed 
to a quasi-isothermal state as seen previously in the 2D 
studies in PS with adiabatic boundary conditions. Thus, 
the source of free energy, the temperature gradient, is 
quickly exhausted in the unstable region before it can be 
replenished by the slow heat conduction from the outer, 
stable layers. We refer to this method of saturation as 
Case A: temperature gradient relaxation. 

An interesting insight into the saturated state can be 
gleaned from examining the average angle of the mag- 
netic field with respect to the horizontal in the unstable 
layer, defined as 

^K sin_1 (pD)- (26) 

This quantity is plotted in Figure 2] for run Fl. The 
average angle of the magnetic field is initially zero — 
consistent with the initial field being perpendicular to 



the temperature gradient. As the MTI grows exponen- 
tially, there is a rapid rise in the angle of the magnetic 
field at the same time as the rapid growth in the kinetic 
and magnetic energy densities. After this initial rise in 
magnetic field angle, (0b) slowly reaches an asymptotic 
value of 18.9° . The saturated state for this case is primar- 
ily reached by temperature gradient relaxation, although 
the effects of the rearrangement of the magnetic field are 
conflated with the former. One can calculate a rough 
steady state heat flux in the anisotropic layer (region 2) 
as, 

Qaniso = — Xaniso.O Sin(^s) [ — — ) 

\ dz J2 

= -(1 x lCT 4 )(sinl8.9 o )(-0.0367) 

« 1.2 x 1(T 6 , (27) 

using the values from Tabled] We can compare the mea- 
sured Coulombic heat flux to a fiducial heat flux. The 
heat flux that would occur if heat were to flow across the 
domain with purely isotropic conductivity for the given 
initial temperature gradient is defined as 

/ AT \ 

Q = -min[xaniso,Xiso] ^— J , (28) 

which is for this run limited by the conductivity of the 
isotropic layer to be Q — 3.3 x 10~ 6 . Thus the steady- 
state anisotropic heat flux is a bit more than one-third of 
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Fig. 1. — Time evolution of the components of the volume-averaged magnetic and kinetic energy density for the entire domain for Run 



Fl, 



1000 3000 3000 4000 
Time (t/rj 



1000 3000 3000 4000 
Time (t/T ) 



Fig. 2. — Evolution of the Alfven velocity and RMS Mach number 
of run Fl averaged over the central MTI-unstable region. The two 
are almost in equipartition. 



the fiducial value. The anisotropic heat flux can be com- 
pared to the steady state heat flux in the isotropic lay- 
ers (regions 1 and 3) which is Qi so = Xiso,o (dT/dz) 1 3 w 
1.3 x 10 -6 . Thus, the steady state heat fluxes are roughly 
consistent. In addition to the Coulombic (eq. [5]) and 
radiative (eq. [7]) heat fluxes, we define an advective 




Fig. 3. — Vertical profile of the horizontally-averaged tempera- 
ture profile in run Fl at various times. Times are labeled in terms 
of t s . The final, saturated state (dotted line) is almost isothermal 
in the central region. 

heat flux in the vertical direction, 

V z e, 



idv 



(29) 

where, V z is the vertical component of the fluid velocity, 
and e = p/("f — 1), is the internal energy. The heat 
fluxes for this run are plotted in Figure [5l The advective 
heat flux is smaller in magnitude than the conductive 
flux and asymptotes to zero from negative values as the 
atmosphere rearranges to an isothermal state. 

The saturation of Case A can be described by a quick 
relaxation of the temperature gradient followed by mod- 
est magnetoconvection until the dynamo has increased 
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Time (t/r ) 

Fig. 4. — Evolution of the average angle of the magnetic field 
with respect to the horizontal in degrees for run Fl. The sharp 
initial peak is followed by a much slower increase to its asymptotic 
value. 

the magnetic field until it is in equipartition with the 
turbulent convective motions. The convective velocities 
are highly subsonic and the magnetic field remains prin- 
cipally confined to the unstable region. More quantita- 
tively, the time for a fluid element to cross the unstable 
layer at the RMS velocity is L/{v) ~ 1.8 x 10 4 r s , which 
is far longer than the time it takes for the simulation 
to saturate. Figure [6] shows a visualization of the mag- 
netic energy density at time, t/rs ~ 3, 840. The buoyant 
plumes are largely confined to the stable region with little 
flux pumping to the stable region. 

Runs Nl, N4, and N5 also saturate with a case A mech- 
anism. N2 is the exact same as the fiducial case except it 
is initialized with a net magnetic flux. The isotropic ther- 
mal conductivity in regions 1 and 3 and the anisotropic 
conductivity in region 2 are a factor of five higher in runs 
N4 and N5 than run Fl, thus facilitating a comparison 
of the role of thermal conductivity Runs N4 and N5 dif- 
fer in that the latter has a net magnetic flux threading 
the domain; whereas, the former has zero net magnetic 
flux. It is interesting to note that in two dimensions, the 
zero net magnetic flux condition requires the magnetic 
energy d ensity to decay away by Cowling's anti-dynamo 
theorem ()Cowling 1119341 ). In three dimensions, magnetic 
dynamo action is permitted, and runs Fl and N4 show 
a sustained increase in the magnetic field. By examining 
Table [5] one can compare the runs at the fiducial time of 
t/r s = 2,050. At this time, the temperature gradient is 
essentially fully relaxed, but the magnetic dynamo has 
not had quite enough time to fully amplify the magnetic 
energy density to equipartition with the kinetic energy. 

Comparing runs Fl and N4, we see that the higher 
conductivity only marginally speeds the increase of the 
magnetic energy density. A slightly larger increase is 
seen in the saturation of the temperature gradient (10.7% 



versus 15.8%). The higher conductivity does result in 
an almost three-fold increase in the RMS Mach number. 
On the whole, the uniform increase in conductivity does 
not qualitatively change the saturation mechanism. The 
saturated heat flux is very similar with {Q)/Q ks 44%. 

Comparing runs Fl versus Nl and N4 versus N5, we 
see that a net magnetic flux threading the domain pro- 
duces a marked increase in the strength of the magnetic 
dynamo action bringing the magnetic energy closer to 
equipartition faster in the runs with the net flux. The 
magnitude of the convective velocities remains the same 
in both cases. The degree of saturation of the temper- 
ature gradient is the same, but the average angle of the 
magnetic field with respect to the horizontal is greater in 
the net flux case. At this early stage of saturation, the 
advective heat flux is often negative as the atmosphere 
reaches the new equilibrium. 

3.2. Case B: Reorientation of Magnetic Field 
Geometry 

We proceed to discuss a very different saturated state 
for the MTI, exemplified by run F3. In this run, the 
conductivity in the isotropic, stable layers is an order of 
magnitude larger than the anisotropic conductivity. Fig- 
ure [7] shows the time evolution of the magnetic energy 
density and kinetic energy density over the entire com- 
putational domain. There is a very strong amplification 
of the magnetic energy density in the unstable domain 
by a factor of approximately 5.3 x 10 3 or 7.8 x 10 3 for 
(B) 2 and (B 2 ) 1 respectively, as defined in Equation [231 
Figure [5] shows the corresponding time evolution of the 
Alfven velocity and the RMS Mach number. The Mach 
number of the turbulence has a brief spike corresponding 
to the initial growth of the instability before reaching a 
very stable 1.8 x 10~ 3 . The Alfven velocity is within a 
factor of two of the turbulent velocity with slight varia- 
tion depending on which average of the magnetic field is 
utilized. 

Up until this point in the discussion, the magnetocon- 
vection in run F3 simply resembles a more energetic ver- 
sion of run Fl. The differences are exposed by exam- 
ing the evolution of the angle of the magnetic field with 
the horizontal, shown in Figure [9l and the evolution of 
the horizontally-averaged, vertical temperature profile, 
shown in Figure 1101 It is obvious that the temperature 
profile evolves quickly to its final value, which deviates 
only slightly from the initial profile, i.e. approximately 
91% of the initial temperature gradient is the asymptotic 
value. In contrast, the average angle of the magnetic field 
evolves quite dramatically to a primarily vertical field 
with (9b) ~ 55°. In this case, which we call Case B, the 
MTI saturates not by temperature gradient relaxation, 
but instead by a reorientation of the magnetic field to a 
primarily vertical field. 

We turn our attention to its transport properties. Fig- 
ure [TT] plots the time evolution of the vertical compo- 
nents of these three heat fluxes and the total heat flux 
at the midplane, and at 80% of the height of the do- 
main in Run F3. At the midplane, the Coulombic heat 
flux is the dominant component of the heat flux with 
a smaller positive contribution from the advective heat 
flux. It is interesting to estimate the efficiency of the 
MTI in transporting heat. We can calculate the compo- 
nents of the heat flux a posteriori from our data. Using 
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Fig. 5. — Time evolution of the horizontally-averaged heat flux at the midplane (left panel) and at 80% height (right panel) of the 
simulation domain in run Fl. The total heat flux (solid line) is subdivided into conductive (dot-dash line), radiative (long dashed line), 
and advective (dotted line) components. 




Fig. 6. — Pseudocolor plot of the magnetic energy density at 
time, t/r s Ri 3,840, in run Fl. Note that the unstable region ex- 
tends from 0.067 < z < 0.13 as denoted by the solid lines. The 
initial magnetic energy density is B 2 = 10~ 10 in the unstable re- 
gion and zero elsewhere. The highly subsonic turbulence drives 
little penetrative convection into the stable region. 



Equation [27] to estimate the heat flux is not appropriate 
here since sin(0B) ^ {sinOs) for large angles. The ideal 
heat flux for comparison purposes is essentially purely 
isotropic transport across the unstable layer, which is 

Q = -Xaniso.o (^f) « (1 x HT 4 ) (0.333) « 3.3xf(T 5 . 

V dz / initial 

(30) 

Thus, the essentially constant Coulombic flux of 1.6 x 
10~ 5 is almost exactly half of the fiducial heat flux. This 
fact leads us to conclude that the MTI can be very effi- 
cient at heat transport. 

One fact that bears further elucidation is the small 
time-average of the advective component of the heat flux. 
With the exception of gravity, our domain is essentially 
symmetric from top to bottom. In addition, the box size 
is small compared to the scale height, so our simualtions 
are intrinsically local in nature. In a real plasma of inter- 
est, the Spitzer conductivity, which varies as n sp ex T 5 / 2 , 
would have a substantial variation in thermal conductiv- 
ity in the vertical direction. Such a conductivity vari- 
ation, combined with a potential variation in optically- 
thin cooling and a global temperature profile, could pro- 
duce a larger net advective flux. Such a realistic config- 
uration will be pursued in future work. 

We now discuss one further aspect of run F3 related to 
the magnetoconvection. Figure [12] shows a visualization 
of the magnetic field energy density in the simulation 
at time, t/r s » 4,700. At this intermediate time, the 
maximum magnetic energy density exceeds B 2 = 10 , 
five orders of magnitude larger than the initial value, 
for regions colored red in the plot. In addition, one 
can see buoyant plumes penetrating into the stable re- 
gion as a result of convective overshoot. Much magnetic 
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Fig. 7. — Time evolution of the components of the volume-averaged magnetic and kinetic energy density for the entire domain for Run 
F3. 
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Fig. 8. — Evolution of the Alfven velocity and RMS Mach number 
of run F3 averaged over the central MTI-unstable region. The two 
are almost in equipartition. 



field has now been stored in the stable layers which ini- 
tially had zero magnetic held (colored blue in the plot). 
The phenomena of penetrative convection and overshoot- 
ing has been thoroughly studied in the magnetocon- 
vection literature of stellar atmospheres (|Tobias et al. I 



Fig. 9. — Evolution of the average angle of the magnetic field in 
the central region with respect to the horizontal in degrees for run 
F3. The sharp initial peak is followed by a steady increase to the 
quite large asymptotic value of (9b) ~ 55° 



I200H: IBrummell. Clune. k Toomre I l2002n . In compar- 
ison to this run, Case A saturation has far less pene- 
trative convection due to the weaker magnetoconvection 
(compare to Figure [6]). 
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Fig. 10. — Vertical profile of the horizontally-averaged tempera- 
ture profile in run F3 at various times. Times are labeled in units 
of t s . The temperature evolution happens almost entirely in the 
first 650 sound crossing times or so. The three times in brackets 
are all essentially collinear profiles. 

Finally, we compare run F3 to the other Case B runs 
(N2, N3, and N6). N2 has the same parameters as F3, 
but has a net magnetic flux. N6 has a fivefold increase in 
both components of the thermal conductivity and a net 
flux. We can compare the fiducial time of t/r s = 2,050, 
which is not yet fully saturated for either simulation. As 
we would expect based on our prior discussion, runs N2 
and N6 have amplified the magnetic field modestly more 
in the same amount of time; however, the temperature 
gradient saturation is similarly minimal in both cases. 
Run N3 is a check for convergence having doubled the 
resolution in each direction. As would be expected for 
a higher resolution simulation, run N3 is approaching 
saturation more quickly than the lower resolution fidu- 
cial case. The asymptotic behavior though results in 
very similar magnetic field amplification, levels of turbu- 
lence, temperature profile saturation, and most impor- 
tantly heat flux. Thus, we have a check that our fiducial 
resolution is sufficient. 

3.3. Comparison of the Cases 

Having outlined the two different saturated states of 
the MTI, we now contrast the two. Table [2] lists the sat- 
urated values of Case A (Run Fl) and Case B (Run F3). 
Case A, temperature gradient relaxation, is quite similar 
in nature to the 2D simulation with adiabatic boundary 
conditions in PS. A small amount of heat flux is supplied 
to the unstable layer, leading to the isothermal temper- 
ature profile. In addition, the field remains primarily 
horizontal with heat flux dominated instantaneously and 
in the time average by thermal conduction. The convec- 
tion is highly subsonic and little magnetic flux pumping 
is observed. The net heat flux is roughly one-half of the 
fiducial value in the saturated state. 

Case B, on the other hand, has a large heat flux sup- 
plied to the unstable layer, which causes the MTI to sat- 
urate by re-orienting the field lines in a largely vertical 
direction. The level of convection, while still subsonic, is 
far more vigorous than Case A. The advective heat flux 
is consistently positive; however, the total heat flux is 



still dominated by the Coulombic component. The con- 
vective plumes are able to traverse the unstable layer in 
approximately 700 sound crossing times, thus leading to 
a significant amount of convective overshoot and mag- 
netic energy storage in the stable layers. The net heat 
flux is roughly one-half of the fiducial value. 

In both cases, the magnetic energy density increases 
through dynamo action until it is roughly in equipartition 
with the kinetic energy density of the subsonic convec- 
tion. It should be noted that the plasma energy density 
is still dominated by the thermal component. In both 
cases, the effective thermal conductivity, as measured by 
the ratio of the actual heat flux to the fiducial heat flux, 
is quite high. 

3.4. The Intermediate Case 

Intermediate to the two cases previously discussed lies 
run F2 where the thermal conductivity in the middle and 
outer layers is constant, Xaniso,o = Xiso,o = 10~ 4 . Figure 
U3l shows the time evolution of the volume-averaged com- 
ponents of the magnetic energy density and kinetic en- 
ergy density over the entire domain, including the stable 
regions. Appreciable amplification of the magnetic en- 
ergy density is seen, especially in the z-componcnt. The 
kinetic energy, after an initial transient, settles down to 
stable constant turbulence with an RMS Mach number 
of 6.1 x 10~ 4 averaged over region 2. The (essentially) 
saturated Alfven velocity is 1.4 x 10~ 4 , thus bringing the 
magnetic energy density into quasi-equipartition with the 
kinetic energy density. The angle of the magnetic field 
with respect to the horizontal saturates at an angle of ap- 
proximately 37°. The heat fluxes for this run are shown 
in Figure [14] The saturated Coulombic heat flux in both 
layers is found to be approximately 9 x 10~ 6 , slightly less 
than one-third of the maximum possible value — still con- 
sistent with a relatively efficient transport of heat. The 
advective heat fluxes, while initially zero, are asymptoti- 
cally oscillating around zero. Since we define the fiducial 
heat flux with respect to the minimum thermal conduc- 
tivity in the problem (imagine a series of thermal resis- 
tors), the intermediate case appears to be less efficient 
than the other two cases. In run Fl, the minimum con- 
ductivity is in the isotropic stable layer, thus causing Q 
to be an order of magnitude smaller for it. If we were 
to instead calculate the fiducial heat flux with only the 
thermal conductivity of the anisotropic layer, the inter- 
mediate case would indeed have a thermal efficiency in- 
termediate to runs Fl and F3. 

Having outlined three different cases of saturation, we 
now show that there is a coherent and logical ordering 
of saturation parameters among them. As we progress 
from Case A to Case B, the free parameter we increase 
is the ratio of isotropic conduction in the outer layers 
to the anisotropic conduction in the unstable region, 
Xiso,o/Xaniso,o- This variation is tantamount to increasing 
the applied heat flux to the MTI-unstable layer. The re- 
sultant saturation properties are evident by simply scan- 
ning Table O Firstly, the vigor of the resultant mag- 
netoconvective dynamo, as characterized by amplifica- 
tion of the magnetic field (6(B 2 ) and 5(B) 2 ), the RMS 
Mach number, and the Alfven velocity, increases signifi- 
cantly as we progress from Case A to Case B. Secondly, 
the resultant atmospheric state changes smoothly in this 
progression as evidenced by the increase in the average 
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Fig. 11. — Time evolution of the horizontally-averaged heat flux at the midplane (left panel) and at 80% height (right panel) of the 
simulation domain in run F3. The total heat flux (solid line) is subdivided into conductive (dot-dash line), radiative (long dashed line), 
and advective (dotted line) components. Both the time-average conductive and advective heat fluxes are positive, but the conductive flux 
dominates. 



magnetic field angle. The degree of relaxation of the tem- 
perature gradient, as shown in Figure [T5l decreases sub- 
stantially from Case A to Case B. Yet, the net heat flux, 




Fig. 12. — Pseudocolor plot of the magnetic energy density at 
time, t/r a m 4,700, in run F3. Note that the unstable region ex- 
tends from 0.067 < z < 0.13 as denoted by the solid lines. The ini- 
tial magnetic energy density is B 2 = 10 — 10 in the unstable region 
and zero elsewhere. Note the plumes due to penetrative convection 
and the storage of magnetic flux in the stable layers. 



which is a convolution of these two properties, monoton- 
ically increases from Case A to Case B. All three cases 
exhibit quite efficient transport of heat flux compared to 
the fiducial value. 

4. APPLICATIONS 
4.1. Plasmas in Clusters of Galaxies 

The most promising potential application of the MTI 
is to the intracluster medium of clusters of galaxies. The 
plasma has a typical density of 10 _3 -10~ 2 cm~ 3 and a 
temperature of roughly 1 — 15 keV, leading to copious X- 
ray emission measured by Chandra and XMM-Newton. 
The total cluster luminosities are found to be typically 
L x ~ 10 43 erg/s to 10 46 erg/s in the X-ray. Most of 
these clusters appear to be in hydrodynamic equilibrium, 
as the time since the last major m erger is approximately 
5 Gyr (jPeterson &: Fabian 1120061 ). although a few excep- 
tional examples have relatively disturbed shapes due to 
recen t major merge rs, such as Perseus (|Fabian et al. I 
[20061) Abell 3667 (|Vikhlinin et al. 1 |2001|), and Abell 
2142 (jMarkevitch et al. 1 120001 ) . Estimates for magnetic 
field strength in clusters range from 1 - 1 /j,G at the cen- 
ter and 0.1 - 7 /iG at a radius of 1 Mpc (jBrunetti et al. I 
12001b iFeretti 1 119991 ) , values which are typically not dy- 
namically important. 

To understand the potential for the MTI, begin by con- 
sidering the Spitzer conductivity for a cluster plasma, 

30 , /fcsT\ 5/2 /10- 2 cm- 3 \ / 37 \ cm 2 

(31) 

where n e is the electron dens ity and In A^ is the Coulomb 
logarithm (|Chandran | [2005"l . Since p e <C A mfp , the per- 
pendicular electron conductivity may be safely ignored, 
e.g. for a 1 fj,G field, k±/ku ~ 10~ 24 . The gravitational 
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Fig. 13. — Time evolution of the components of the volume-averaged magnetic and kinetic energy density for the entire domain for Run 
F2. 
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Fig. 14. — Time evolution of the horizontally-averaged heat flux at the midplanc (left panel) and at 80% height (right panel) of the 
simulation domain in run F2. The total heat flux (solid line) is subdivided into conductive (dot-dash line), radiative (long dashed line), 
and advective (dotted line) components. The conductive flux clearly dominates the advective flux. 
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Fig. 15. — Comparison of the horizontally-averaged vertical tem- 
perature profiles. The initial condition is plotted with a solid line. 
The saturated temperature profiles of run F3 (dotted line), run F2 
(dashed line), and run Fl (dash-dot line) show increasing degrees 
of relaxation of the temperature profile towards isothermality in 
the unstable layer. 

potential of the cluster is well -approximated by an NFW 
profil e for the dark matter (jNavarro. Frenk. fc White! 
Il997h . Now we can consider the instability criterion for 
the MTI, Equation [T21 for a typical cluster with mass 
10 14 M , an NFW scale radius of 500 kpc, and a typical 
density of 10~ 2 g/cm 3 . The first term of the instability 
criterion is a measure of magnetic tension at a given scale 
and field strength. For our model cluster, the plasma 
beta is (3 s=s 2, 000, thus, the magnetic field is dynamically 
weak. The second term is set by the product of the clus- 
ter gravitational cluster potential (p~ 1 (dp/dz) = V<E>) 
and the logarithmic derivative of the temperature. We 
calculate a reasonable temperature profile by using a 
fit to cluster struct ure formation calculations given by 
lLoken et al.~l (|2002f) . Finally, we assume that outside of 
the cluster core, conduction is entirely dominated by par- 
allel conduction, thus x'c/x' = 1- Balancing these two 
terms gives the smallest wavelength, A cr it, at which the 
ICM is marginally unstable to the MTI. This lengthscale 
can be parameterized as 



A„it = 4.6 kpc 



T 



5 keV 



1/2 



/2000\ 



1/2 



(32) 



Thus all scales larger than 4.6 kpc are unstable to the 
MTI. An alternative way of looking at this problem is 
to estimate what field strength would be required for 
the largest scale, i.e. the cluster virial radius, to be 
stable to the MTI. For our model cluster this requires 
« 0.16, which corresponds to magnetic field strengths 
in excess of 110 /xG. Typical measurements of cluster 



field strengths by the various mechanism of synchrotron 
emission, Faraday rotation, or in verse Compton emissio n 
never approach this large value (|Carilli fc Taylor 1 12002') . 

In addition, we can comment on the temperature pro- 
files of clusters. The previously mentioned structure for- 
mation calculations assume a ACDM cosmology and pre- 
dict a moderately steep, monotonically de creasing tem- 
perat ure profile of the intracluster plasma (|Loken et al. I 
2002). Recent observations with XMM-Newton have in- 
ferred a temperature profile that is flatter than the pre - 
dictions from structure formation (|Pratt et al. I [2006) . 
It may be that the observed temperature profiles are in- 
dicative of the saturated state of the MTI. The MTI also 
provides for a dynamo mechanism to amplify a small seed 
magnetic field through magnetoconvection. 

A number of authors have considered the effect 
of thermal conduction on the temperature profile 
of clusters of galaxies using an effective conductiv- 
ity across tangled field l i nes in a turbulent medium 
(iRechester fc Rosenbluth I 119781 : IChandran fc Cowlevl 
1998). More recent efforts have revised this calculation 
to include a spectrum of correlation lengths, resulting in 
an effective conductivity that is much closer to the value 
of the full Spitzer conductivity dNaravan fc Medvedev I 
120011: IZakamska fc Naravan 112003^ " In this case the tur- 
bulence which tangles field lines is externally driven. 
However, the ICM is unstable to the MTI, and therefore 
the turbulence, the tangling of the field lines, and the re- 
sulting heat flux are determined by the same underlying 
driving mechanism. Computing the temperature profile 
and saturated heat flux self-consistently in a typical ICM 
is an important goa l of future efforts. 

iChandranl (|2005h has constructed an analytical model 
that includes the effects of anisotropic cosmic ray pres- 
sure in addition to the thermal conduction in the MTI 
model. It is worth noting that non-thermal particles 
may also contribute t o the structure and dynam ics of 
the ICM. In addition, IChandran fc Dennis I ((20061 ) have 
proposed that a cosmic ray version of the MTI could be 
at play in the core of the cluster, where a radial gradient 
of cosmic ray pressure (from an AGN for example) can 
play the role the temperature gradient does in the MTI. 
This heat transport in the core is a potential mechanism 
to obviate the need for cooling flows. 

4.2. Other Applications 

Two other potential applications of the magnetother- 
mal instability have also been identified. First, neutron 
star atmospheres have large temperature and pressure 
gradients due to a surface gravity of g « 2 x 10 14 cm/s 2 , 
giving a scale height of around 10 cm. The temperature 
gradient must be large to drive the MTI in the presence 
of large magnetic fields. For example, millisecond radio 
pulsar s typically have magn etic fields in the 10 8 - 10 9 G 
range (jBhattacharva 1 119951 ). In a region near the equa- 
tor and inside of the photosphere, there is the potential 
for the MTI to be unstable. Here the primary limiting 
factor is the term x'c/x' from the instability criterion. 
The atmosphere is a semi-degenerate, semi-relativistic 
plasma. Deep in the atmosphere, the Larmor radius be- 
comes comparable to the mean free path. Near the sur- 
face of the atmosphere, radiative transport primarily due 
to Brehmsstrahlung becomes important. Both of these 
properties suppress the MTI due to isotropic transport. 
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Between these two regions, though, the neutron star at- 
mosphere may be MTI-unstablc. Due to the complex 
nature of the physics, this supposition will require de- 
tailed modeling of the neutron star atmosphere to assess 
with certainty. 

A second potential application is to Radiatively Inef- 
ficient Accretion Flows (RIAFs). Observations of the 
electron mean free path in Sgr A* in our Galactic Center 
show that A m f p ~ 10 5 - 10 6 Schwarschild radii. Thus, 
the plas ma is strongly collisionless. Recent analytical 
work bv lJohnson &; Quataert I (|2006h shows that includ- 
ing electron thermal conduction can suppress the accre- 
tion rate by 1 - 3 orders of magnitude from the Bondi 
rate — a first step in explaining the very low luminosity of 
Sgr A*. This analysis is static, however, and cannot fully 
account for the resul tant buoyant m o tion d riven by the 
thermal conduction. ISharma et al. I (|2003f ) have devel- 
oped Landau fluid closures that describe kinetic MHD 
appropriate to this regime. Fully collisionless shea ring 
box simulations performed bv ISharma et ai. I ipooeh ac- 
tually lead to a modestly increased a-parameter com- 
pared to typical MRI turbulence in an MHD simulation. 
This increase in a is due to anisotropic pressure effects in 
the fully collisionless treatment. Such an increase is con- 
sistent w ith the magne toviscous instability in the shear- 
ing box (jBalbus Il2004f k however, it cannot account for 
observations of much lower accretion rates for Sgr. A*. 
As a result, global simulations of collisionless accretion 
are called for that include the temperature stratification 
of the flow. 

5. SUMMARY AND FUTURE WORK 

Many astrophysical plasmas fall into a low density, 
weakly magnetized regime, in which transport of heat 
and momentum due to thermal conduction and viscos- 
ity takes place primarily along magnetic field lines. One 
would expect an atmosphere with dS/dz > to be com- 
pletely stable to convection. Yet, the magnetothermal 
instability destabilizes this atmosphere if dT/dz < 0. We 
have performed fully three-dimensional time-dependent 
MHD simulations of the nonlinear regime of this insta- 
bility. 

The nonlinear evolution of the MTI in 3D with differ- 
ent parameters shows several common features. Expo- 
nential growth in the linear growth rate produces a spike 
in all the components of the kinetic and magnetic energy 
densities, and a rapid rearrangement of the temperature 
profile to the final saturated temperature gradient. Af- 
ter this initial transient, the atmosphere enters a mag- 
netoconvective phase in which the fluid motions reach a 
steady level of subsonic turbulence. Simultaneously, the 
magnetic field is amplified in a dynamo until the mag- 
netic energy density is in approximate equipartition with 
the kinetic energy density. 

As was outlined in fj3j we have identified two different 
saturated states for the MTI in 3D. These mechanisms 
differ proximately in the strength of the isotropic thermal 
conduction in the stable layers relative to the strength of 
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the anisotropic thermal conduction in the MTI-unstable 
region. This adjustable parameter allows an exploration 
of the parameter space for a range of driven heat flux 
boundary conditions. In Case A: temperature gradient 
relaxation, the isotropic heat flux is small and the un- 
stable layer saturates in a way reminiscent of adiabatic 
simulations from 2D. In the saturated state, the tem- 
perature profile is almost isothermal and the magnetic 
field is still primarily horizontal in orientation. In Case 
B, the isotropic heat flux is large and the unstable layer 
saturates in a significantly different way. The tempera- 
ture gradient relaxes only slightly; however, the magnetic 
field rearranges itself into a primarily vertical configura- 
tion to carry the applied heat load. This variation of 
the isotropic heat flux changes the vigor of the magneto- 
convection, resulting in higher RMS Mach numbers and 
higher magnetic field amplification as we approach Case 
B. In all cases, though, we find the instability reaches a 
natural quasi-equipartition between subsonic turbulent 
motions and magnetic energy density. 

A key conclusion from this work is that the magnetic 
field can no longer be thought of a passive contami- 
nant. We began with an atmosphere that was convec- 
tively stable to the Schwarschild criterion with a weak 
magnetic field. We can quantify "weak" in terms of the 
plasma beta, ft = Airp/B 2 , the ratio of thermal to mag- 
netic pressure. For our initial choice of magnetic field 
in the simulations, ft — 10 10 . Thus, conventional wis- 
dom leads to considering the magnetic field as dynam- 
ically irrelevant — a passive contaminant simply frozen 
into the background fluid and transported by purely hy- 
drodynamic motions. Yet, the illustrative example of 
Case B shows that this reasoning is patently incorrect 
in many astrophysical plasmas. The transport of heat 
along the magnetic field triggers the magnetothermal in- 
stability which is capable of reordering and amplifying 
the magnetic field, generating magnetoconvection, and 
transporting heat. The dynamically weak magnetic field 
now completely dominates the resulting dynamics. 

The most immediate application for future study is 
to the intracluster medium of clusters of galaxies. A 
fully three-dimensional MHD simualtion including real- 
istic Spitzer thermal conductivity and radiative cooling 
is required. In addition, the effects of conduction on ra- 
diatively inefficient accretion flows and convection domi- 
nated accretion flows must be investigated in the context 
of MRI turbulence in global MHD simulations. 
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